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Abstract 



The Navier-Stokes transport coefficients of a granular gas are obtained 
from the Chapman-Enskog solution to the Boltzmann equation. The granu- 
lar gas is heated by the action of an external driving force (thermostat) which 
does work to compensate for the collisional loss of energy. Two types of ther- 
mostats are considered: (a) a deterministic force proportional to the particle 
velocity (Gaussian thermostat), and (b) a random external force (stochastic 
thermostat). As happens in the free cooling case, the transport coefficients 
are determined from linear integral equations which can be approximately 
solved by means of a Sonine polynomial expansion. In the leading order, we 
get those coefficients as explicit functions of the restitution coefficient a. The 
results are compared with those obtained in the free cooling case, indicating 
that the above thermostat forces do not play a neutral role in the transport. 
The kinetic theory results are also compared with those obtained from Monte 
Carlo simulations of the Boltzmann equation for the shear viscosity. The com- 
parison shows an excellent agreement between theory and simulation over a 
wide range of values of the restitution coefficient. Finally, the expressions of 
the transport coefficients for a gas of inelastic hard spheres are extended to 
the revised Enskog theory for a description at higher densities. 
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I. INTRODUCTION 



The usefulness of fluid-like type equations to describe systems of granular particles in 
rapid, dilute flow has been recognized for many years. The essential difference from ordinary 
fluids is the absence of energy conservation, yielding subtle modifications of the conventional 
Navier-Stokes equations for states with small spatial gradients of the hydrodynamic fields. 
Although many efforts have been made in the past few years in the understanding of these 
systems, the analysis of the influence of dissipation on the transport coefficients still remains 
a topic of interest and controversy. For a low density gas, these coefficients may be deter- 
mined from the Boltzmann equation modified to account for inelastic binary collisions. The 
idea is to extend the Chapman-Enskog method [T[ to the inelastic case by expanding around 
the local version of the homogeneous cooling state (HCS), i.e., a reference state in which 
all the time dependence occurs through the granular temperature. In the first-order of the 
expansion, explicit expressions for the transport coefficients as functions of the restitution 
coefficient have been obtained in the case of hard spheres as well as for a rf-dimensional 
system ||. This analysis has been also extended to higher densities in the context of the 
Enskog equation H. The results obtained in the latter theory describe very well the hydro- 
dynamic profiles obtained in a recent experimental study of a three-dimensional system of 
mustard seeds fluidized by vertical vibrations of the container [[J. 

One of the main difficulties in obtaining the above transport coefficients lies in the fact 
that, in contrast to what happens for molecular fluids, the reference state (zeroth-order so- 
lution of the Chapman-Enskog expansion) depends on time due to the dissipation of energy 
through collisions. In addition, it is also well-known that this state is unstable to long enough 
wavelength perturbations so that the state becomes inhomogeneous for long times ||. To 
overcome such difficulties, one possibility is to introduce external forces to accelerate the 
particles and hence compensate for collisional cooling. As a consequence, the corresponding 
reference state is stationary and linearly stable against spatial inhomogeneities. This mecha- 
nism of energy input (different from those in shear flows or flows through vertical pipes) has 
been used by many authors in the past years to analyze different problems, such as 

non-Gaussian properties (cumulants, high energy tails) of the velocity distribution function 
P]T2||, long-range correlations ||10|| , and collisional statistics and short scale structure []TT 



Since the latter requires the solution of the corresponding linearized hydrodynamic equa- 
tions around the homogeneous state, the explicit expressions for the transport coefficients 
are needed. Given that the dependence of these coefficients on the restitution coefficient is 
not known, the expressions of these coefficients are usually assumed to be the same as those 
for the elastic gas. However, according to the results derived in the free cooling case , 
the above assumption could be only justified in the small inelasticity limit. 

The goal of this paper is to determine the transport coefficients of a heated granular gas. 
This allows us to measure the new effects induced by the external force on transport by 
comparison with the results derived in the unforced case In addition, one could also 

assess to what extent the previous results on short and large structure |]TT],[n]] are indicative of 
what happens for finite degree of dissipation. There are different mechanisms to inject energy 
to the gas. Here, the fluidization is driven by the action of external forces (thermostats) 
acting locally on each particle. In this paper we will consider two types of thermostats: the 
Gaussian and the stochastic thermostats. In the case of the Gaussian thermostat, the gas is 
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heated by the action of an external force proportional to the peculiar velocity. This type of 
"anti-drag" force can be justified by Gauss's principle of least constraints 0] and has been 
widely used in nonequilibrium molecular dynamics simulations of molecular fluids. Another 
mechanism for thermostatting the system is to assume that the particles are subjected to a 
random external force, which gives frequent kicks to each particle between collisions. If this 
stochastic force has the properties of a white noise, it gives rise to a Fokker-Planck diffusion 
term in the Boltzmann equation M. 

It must be remarked that in most experiments, energy is added to the granular gas 
through a boundary, causing gradients in the energy perpendicular to that boundary. In 
this sense, we do not claim that the above forcing terms are the most suited to model any 
particular real system. However, they have the advantage of that it can be incorporated 
into the kinetic theory very easily. This allows for instance, to test the assumptions of the 
Chapman- Enskog method through a direct comparison with computer simulations, such as 
Monte Carlo or molecular dynamics simulations |T^J] . 

The plan of the paper is as follows. In Sec. we review the Boltzmann equation and 
associated macroscopic conservation laws in the presence of the external forces discussed 
above. The Chapman- Enskog method for solving this equation is presented in Sec. |T| 
and subsequently applied to the cases of Gaussian and stochastic thermostats. As happens 
in the free cooling case 0, the transport coefficients are determined from linear integral 
equations which can be approximately solved by means of a Sonine polynomial expansion. 
In the leading order, we get the transport coefficients as explicit functions of the restitution 
coefficient. Section |l| ends with a comparison between the results derived here in both 
driven cases with those previously obtained in the unforced case 0. Such a comparison 
shows that, in general, the thermostats do not play a neutral role in the transport since they 
clearly affect the dependence of the Navier-Stokes transport coefficients on the dissipation. 
In order to check the degree of reliability of the Sonine approximation, a comparison with 
direct Monte Carlo simulation of the Boltzmann equation is carried out in Sec. [TV} More 
specifically, the simulations are performed for a gas undergoing uniform shear flow, using the 
Gaussian and the stochastic thermostats to control inelastic cooling. In the long time limit, 
a (reduced) shear viscosity can be measured in both simulations. The comparison with the 
Chapman-Enskog solution shows an excellent agreement, indicating that the Sonine results 
have an accuracy comparable to that for elastic collisions. The paper is closed in Sec. [V] 
with a brief summary and discussion of the results presented. In addition, the corresponding 
expressions of the transport coefficients for a granular gas of hard spheres in the framework 
of Enskog kinetic theory are also displayed in Appendix [FJ. 

II. HEATED GRANULAR GASES 

We consider a granular gas composed by smooth inelastic disks (d = 2) or spheres (d = 3) 
of mass m and diameter a. The inelasticity of collisions among all pairs is characterized by 
a constant restitution coefficient a < 1. In the low-density regime, the evolution of the one- 
particle velocity distribution function /(r,v;t) is given by the Boltzmann kinetic equation 

E3H§ 

(d t + Vf V + ^)/(r jVl ,t) = J[ Vl |/(t),/(t)] , (1) 
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where the Boltzmann collision operator J[vx\f, /] is 



J[vi\f,f] =o d - 1 Jdv 2 Jd*e(&-g)(a- g ) 

x a~ 2 f(r, v[,t)f(r, v 2 , t) - /(r, v 1; t)f(r, v 2 , t) 



(2) 



In Eq. (|T]), T is an operator representing the effect of an external forcing which injects energy 
into the granular gas allowing it to reach a steady state. Furthermore, d is the dimensionality 
of the system, <x is a unit vector along their line of centers, is the Heaviside step function, 
and g = v x — v 2 . The primes on the velocities denote the initial values {v^v^} that lead 
to {vi, v 2 } following a binary collision: 



Vl 



a 



v 2 



a 



[a- ■ g)(r 



(3) 



The macroscopic balance equations for density n, momentum mu, and energy |nT follow 
directly from Eq. ([!]) by multiplying with 1, mvi, and and integrating over vi: 



D t n + nV • u = , 



(4) 



D tUi + (mn^VjPij = , (5) 

D t T + A (V • q + PtjVjUi) = -(C - A)T . (6) 
In the above equations, D t = d t + u • V is the material derivative, 

P = J dvmVV/(v) (7) 

is the pressure tensor, 

q = y" dv\mV 2 V /(v) (8) 

is the total heat flux, and V = v — u is the peculiar velocity. In the right hand side of 
the temperature equation (|6|), the cooling rate ( (measuring the rate of energy loss due to 
dissipation) and the source term A (measuring the rate of heating due to the external force) 
are given by 



A = -^/rfv ml , 2 F/(v). (10) 

It is assumed that the external driving does not change the number of particles or the 
momentum, i.e., 

J dvJF/(v) = J dvvf/(v) = 0. (11) 
4 



In the case of elastic particles (a = 1) and in the absence of external forcing (JF = 0), it is 
well known that the long-time uniform solution of Eq. (|l|) is the Maxwell-Boltzmann distribu- 
tion function. However, if the particles collide inelastically (a < 1) and T = 0, a steady state 
is not possible in uniform situations since the temperature decreases monotonically in time. 
In this case, Goldshtein and Shapiro [I5[] showed that Eq. (|1|) admits an isotropic solution, de- 
scribing the homogeneous cooling state (HCS), in which all the time dependence of / occurs 

Tm- d / 2 v d {t)<S>{v/v {t)). 



only through the thermal velocity Vo(t) = y2T(t)/m: /(v, t) 
So far, the exact form of <£> has not been found, although a good approximation for thermal 
velocities can be obtained from an expansion in Sonine polynomials. In the leading order, 
$ is given by 




* i 



(d + 2)v* 2 + 



d(d + 2) 



v/uo, 



where the estimated value of c is 



32(1 - a)(l - 2a 2 



9 + 24d - a(41 - 8d) + 30(1 - a)a 2 ' 



(12) 



(13) 



The estimate ( |i~3"D presents a quite good agreement with Monte Carlo simulations of the 
Boltzmann equation 1)121,17 . 



However, by driving a granular gas by boundaries or external fields it can reach a steady 
state. The energy injected in the gas may exactly compensate for the energy dissipated 
by collisional cooling. The same effect can be obtained by means of external forces acting 
locally on each particle. These forces, which we will call thermostats, are represented by the 
operator T in Eq. (|) and depend on the state of the system. Several types of thermostats 
can be used. Here, we will consider two. One of them is a deterministic thermostat widely 
used in nonequilibrium molecular dynamics simulations of elastic particles which is 
based on Gauss's principle of least constraints. In this case, T is given by [12|] 



*7(v) = \^ • [v/(v)] 



(14) 



where, according to Eqs. (§) and (|H]), the thermostat has been adjusted to get a constant 
temperature in the long time limit. It must be pointed out that the corresponding Boltzmann 
equation (|l|) for this Gaussian thermostat force is formally identical with the Boltzmann 
equation in the HCS (i.e. with T = 0) when both equations are written in terms of the 
reduced distribution $(t>*). As a consequence, the result (|l|) applies to this thermostatted 
case as well. 

Another way of heating the gas is by means of a stochastic force assumed to have the 
form of a Gaussian white noise |7|. The corresponding operator T has a Fokker-Planck form 



(15) 



where again the strength of the correlation has been chosen to achieve a time independent 
temperature. By using this thermostat, van Noije and Ernst have studied the stationary 
solution to the uniform Boltzmann equation ([]]) and found for the coefficient c the value 
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32(l-«)(l-2a 2 ) 



73 + 56^ - 3a(35 + 8d) + 30(1 - a)a 2 



As happened in the Gaussian case, recent Monte Carlo simulations of the Boltzmann equa- 
tion |12|] agree quite well with the Sonine estimate ([II 



III. CHAPMAN-ENSKOG SOLUTION: NAVIER-STOKES TRANSPORT 

COEFFICIENTS 

As said in the Introduction, our goal is to get the Navier-Stokes transport coefficients 
in the presence of the thermostats introduced in the previous Section. This allows one 
to assess the influence of these thermostats on transport by comparison with the results 
obtained in the free cooling case [@||. To do that, we consider now a spatially inhomogeneous 
state, created by initial preparation or by boundary conditions. We assume that the spatial 
variations of n, u, and T are small on the scale of the mean free path. Under these conditions, 
the Chapman-Enskog method |l[ provides a solution to the Boltzmann equation based on 
an expansion around the local version of the heated homogeneous state induced by the 
thermostat forces. This is obtained from the above states by replacing the temperature, 
density, and flow velocity by their nonequilibrium local values. As a consequence, the local 
version of the operators T consists of replacing v ^ V = v — uin Eqs. ([14]) and (p~5|) . 

The Chapman-Enskog method assumes the existence of a normal solution in which all 
the space and time dependence of the distribution function appears through a functional 
dependence on the hydrodynamic fields 

f(v,v,t) = f[v\n(t),u(t),T(t)}. (17) 

For small spatial variations, this functional dependence can be made local in space and time 
through an expansion in gradients of the fields. To generate the expansion, it is convenient 
to write / as a series expansion in a formal parameter e measuring the nonuniformity of the 

system, 

/ = / (0) +e/ (l) +e 2 /( 2) + ... ) (lg) 

where each factor of e means an implicit gradient of a hydrodynamic field. The local reference 
state is chosen such that it has the same first moments as the exact distribution /, or 
equivalently, the remainder of the expansion must obey the orthogonality conditions 

/dv[/(v)-/W(v)]=0, /dw[/(v)-/^(v)]=0, (19) 

^ 2 [/(v)-/ (0) (v)] =0. (20) 

The time derivatives of the fields are also expanded as d% = 9 t + ed^ + ■ ■ •. The coef- 
ficients of the time derivative expansion are identified from the balance equations (H)-(§|) 
with a representation of the fluxes, the cooling rate £ and the heating term A in the macro- 
scopic balance equations as a similar series through their definitions as functionals of the 
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distribution /. This is the usual Chapman- Enskog method for solving kinetic equations [|TJ. 
The main difference with respect to the unforced case @ is that now (as happens in the 
elastic case) the sink term in the energy equation is zero, so that the terms coming from the 
time derivative 9 t vanish. 

Now, we derive the corresponding hydrodynamic equations in the presence of the Gaus- 
sian and the stochastic thermostats. 



A. Gaussian thermostat 

In the case of the Gaussian thermostat fll4]) , the Boltzmann equation becomes 



d t f + v • V/ + • (V/) = J [f, f] . (21) 

Substitution of the Chapman-Enskog solution ([18]) into Eq. fl2~T|) leads to different kinetic 
equations for the distributions f( k \ To zeroth order in e, the Boltzmann equation reads 

k (0) 4;-(Vf {0) ) = j\f {0) ,f {0) ], ( 22 ) 



2^ <9V 

where use has been made of the macroscopic balance equations at this order 

d^n = 0, d? } u = 0, «9 t (0) T = 0. (23) 
Here, the cooling rate is determined by Eq. (|9|) to zeroth order 

C (0) = -^ / dvmv 2 J[f(°\f(% (24) 

The solution to Eq. (|22"D /( ) = f^°\V) is isotropic so that the zeroth order approximations 
to the pressure tensor and heat flux are 

P§*=p8a, q (0) = 0, (25) 

where p = nT is the hydrostatic pressure. The distribution /(°) is essentially given by the 
Sonine approximation QT2] ) with the cumulant c given by Eq. flT3|). 

The analysis to first order in e is similar to the one woked out in Ref. for the free 
cooling case. We only display here the final expressions for the fluxes with some details 
being given in Appendix 0. The final result to first order in the spatial gradients is 

PiP = -V hiUj + VjUi - -dijV ■ u) , (26) 



q« = - K VT - fiVn, (27) 

where rj is the shear viscosity, k is the thermal conductivity, and /i is an additional transport 
coefficient not present in the elastic case. In dimensionless form, the transport coefficients 
are given by 
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(28) 



K 



K 
K 



d-1 1 + c 



d 



(29) 



Tkq 



d-1 
2d v* 



A* 2^ 



(30) 



Here, r] and k are the elastic values of the shear viscosity and thermal conductivity [|TJ 
respectively, 



Vo = ^n-^T (d/2) (mT)^a-^\ n = 

8 2{d — 1) m 



(31) 



C* = C /^o, = p/^o is a characteristic collision frequency, and c(a) is related to the 
fourth moment of by 



c(a) 



8 



d(d + 2) 



277 nJ J 



d{d+2) 



(32) 



The coefficient c(a) is a measure of the deviation of the reference state from that of a gas with 
elastic collisions. As said above, a good estimate of c(a) is given by the Sonine approximation 

, we have introduced the reduced collision frequencies 



(|l3l). Furthermore, in Eqs. (|2£ 



//,, 



i/o/<*vAi(V)Cy(V) : 



(33) 



JrfvS(V)-£A(V) 
' K i/„/dvS(V)-.A(V) 



/rfvS(V) ■ £g(V) 
zo/dvS(V)-8(V)- 



where 



A,(V) = m ( VSV, - ^V 2 ) , S(V) = { > 2 - d -±^ 



~2~ 



T V. 



(34) 



(35) 



So far, all the results are exact but not explicit because of £*, i/* i/*, and To get more 
explicit expressions for the dependence of the transport coefficients on a it is convenient to 
use the leading Sonine polynomial approximations for A(V), B(V), Cij(V), and f^°\ In the 
case of /(°', we take the approximation (|j"2"|) while the remaining quantities are given by 



/ A(y) \ 




( ctS(V) \ 




B(V) 


- f M (V) 


CnS(V) , 


fu{V) = nix 






vCn A,(v); 





-(V/v ) 2 



(36) 



The factor fm(V) occurs since these polynomials are defined relative to a Gaussian scalar 
product. The coefficients are the projections of A, B, and Cij along S(V), and -Dy(V), 
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C T 



ni 



d(d + 2) nT 3 



{ B(V) • S(V) 



' d+2 nT 2 ^ 
__2_rn_ n 
d+2 T 3 r 



(37) 



1 



(d + 2)(d 



T)^F 2 I dYC ^) D ^ v ) 



nT 2 



T) . 



(38) 



where use has been made of the definitions ( |A11| ), ( |A12| ), and ( |A13|) . With these expressions, 
the cooling rate (* and the collision frequencies v*, i/*, and v* can be explicitly computed. 
These calculations have been made by Brey and Cubero 0, with the result 



d + 2 



a 



32 



(39) 



(l -a + |d) (1 + a) (l - 



c 

64 



(40) 



v.. 



1 + a 



d-1 3 . , oWl , 4 + 5d-3(4-d)a 

1 (d + 8)(l - a) H — c 

2 16 V A ; 1024 



(41) 



where c is given by Eq. 

Substitution of the above expressions into Eqs. (|28|)-(|30D gives finally the explicit de- 
pendence of the transport coefficients on the restitution coefficient a. In order to gain some 
insight into the behavior of 7]*, k*, and fi* it is convenient to consider the weak dissipation 
limit. In this limit, the above coefficients can be expanded in powers of the inelasticity 
parameter 1 — a 2 . The leading contributions are 



QAd 2 - 97d + 32 
128d(d-l) ' 



« 2 ) + 



(42) 



1 + 



56d- 191 
128(d- 1)' 



a 2 ) + 



(43) 



4(1 -d) 



(1 



a 2 ) + 



(44) 



B. Stochastic thermostat 

In the case of the stochastic thermostat ([IB]) , the Boltzmann equation becomes 

9t/ + v .V/-^ C (^) 2 / = ^[/,/]. (45) 

As before, this equation can be solved by means of the Chapman-Enskog method. Since the 
procedures to arrive at the expressions of the transport coefficients are identical to the ones 
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employed in the case of the Gaussian force, here we only quote the final expressions of the 
transport coefficients. In dimensionless form, they are 



(46) 



d- 1 1 + c 
d v* 



(47) 



d — 1 c 

2d ^* 



(48) 



Here, in the first Sonine approximation, the expressions of (*, u* and u* are also given by 
Eqs. ([39p-pT|), respectively, with c(a) given by Eq. (0). In the quasielastic limit, one has 



128d(rf-l) {1 ~ a) 



(49) 



„ 1 8d- 287 ^ 2 , 

+ i~28(rf-^T) ~ a ) + ' 



(50) 



4(1 -d) 



;i-a 2 ) + - 



(51) 



C. Comparison with the free cooling case 

The transport coefficients for a dilute unforced granular gas in d dimensions has been 
recently obtained by Brey and Cubero ||. These authors generalize a previous derivation 
made from the Boltzmann equation in the three dimensional case The results are 



(52) 



d - 1 1 + c 

d u* - 2C* ' 



(53) 



2(* 



2< - 3C* 



k + 



(d- l)c 
2d(* 



(54) 



where c(a) is given by Eq. ([13]). In the quasielastic limit, these coefficients behave as 

, 1 48^ -113d + 64.., 2 . 
*7 -1 + — , noJfJ — - — (l-a 2 ) + ---, (55) 



128d(d- 1) 
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(56) 



^s^rV-^- < 57) 



Comparison between Eqs. fl42|)-([44D, (|49|)-(|5T1), and (pB])-(|5Y|) shows that the leading 
order corrections to the elastic values of the transport coefficients are clearly different in 
the free cooling and driven cases. This illustrates the fact that the transport properties are 
affected by the thermostat introduced so that the latter does not play a neutral role in the 
problem. In Figs. [I], |2|, and |3| we plot the reduced coefficients 77*, k*, and //*, respectively, 
for hard spheres (d = 3) in the cases of the unforced gas, the Gaussian thermostat, and the 
stochastic thermostat. In general, we observe that the predictions obtained in the Gaussian 
case for the shear viscosity and the thermal conductivity are closer to those of the unforced 
case than to the ones obtained with the stochastic force. This is in part motivated by the 
fact that, in dimensionless form, the results of the unforced and Gaussian cases are identical 
in the uniform problem (reference state). In the case of the shear viscosity, all the theories 
give the same trends for rj* since this coefficient increases with the dissipation. However, 
at a qualitative level, the influence of dissipation on the viscosity in the stochastic case is 
much less significant than in the other two cases. Thus, for instance, for a = 0.8 (moderate 
dissipation), the shear viscosity of the granular gas rj has only changed about 1% with respect 
to its elastic value when the gas is heated by means of the stochastic force. This could justify 
the use of the elastic values of the transport coefficients in the analysis of the long-range 



correlations in a randomly driven granular fluid [ID||. Discrepancies between unforced and 
the two driven systems are more important in the case of the thermal conductivity k*, as 
Fig. ^| shows. In addition, while in the free cooling and Gaussian cases k* increases with 
a, the opposite happens in the stochastic case. Finally, Fig. |3| shows the dependence of the 
coefficient /i* on a. We observe that /i* ~ for both thermostats, while this coefficient 
is clearly different from zero in the free cooling problem. This is basically due to the fact 
that /i* oc c(a) in the Gaussian and stochastic cases, Eqs. ( |3~0"D and (0), and so it vanishes 
exactly if one takes the Maxwellian approximation (which is known to give a very accurate 
description) for the reference state. This means that, for practical purposes, one can neglect 
the contribution to the heat flux coming from the term proportional to the density gradient 
in the heated gas case: q^ 1 ' — > — kVT. As is apparent from Fig. [| this cannot be assumed 
in the unforced case since /1* is clearly different from zero even in the quasielastic limit. 



IV. COMPARISON WITH MONTE CARLO SIMULATIONS 

As has been discussed above, the practical evaluation of the transport coefficients re- 
quires the truncation of an expansion for the solutions of the integral equations in Sonine 
polynomials. In the case of elastic collisions, the leading order truncation is known to be 
a very good approximation M. A natural question is whether the above degree of accu- 
racy is also maintained in the inelastic case. To answer this question, one has to resort 
to numerical solutions of the Boltzmann equation, such as those obtained from the Direct 
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Simulation Monte Carlo (DSMC) method fl9 |. Although this method was originally devised 



for molecular fluids, its extension to deal with inelastic collisions is straightforward [12,|l8| 

Recently, the shear viscosity of a low-density granular gas of freely evolving hard spheres 
has been determined from the DSMC method [^UJ. The results show a very good agreement 



with the predictions based on the Boltzmann equation in the first Sonine approximation, 
Eq. (p2|). This experiment consists in preparing an initial inhomogeneous nonequilibrium 
state corresponding to a transverse shear wave, and then analyzing its subsequent evolution 
in time. The shear wave decays exponentially with a time scale inversely proportional to the 
viscosity An alternative route to measuring the shear viscosity consists of preparing a state 
of uniform shear flow using Lees-Edwards boundary conditions [^TJ. Macroscopically, this 
state is characterized by constant density n, a uniform temperature T, and a linear velocity 
profile a = du x /dy = const. In a molecular fluid, unless a termostatting force is introduced, 
the temperature increases in time due to viscous heating. The corresponding energy balance 
equation can be used to determine the shear viscosity for sufficiently long times In a 
granular fluid, the relationship between the temperature and the shear viscosity is not simple 



since there is a competition between viscous heating and collisional cooling . However, 

if external forces of the form (|14|) and (|15|) that exactly compensate for the collisional energy 
loss are introduced, the viscous heating effect is still able to heat the system. Under these 
conditions, the Boltzmann equation to be solved is 

dtf ~ aV y ^-f + ^ ■ (V/) = J [f, f] (58) 
in the case of the Gaussian thermostat, and 

d ' ! - aV *wJ-\lc{^)' s-nu] (59) 

in the case of the stochastic thermostat. 

Equations (|58| ) and (|59| ) have been numerically solved by means of the DSMC method for 
a three-dimensional system (d — 3). At given values of the shear rate a and the restitution 
coefficient a, we start from a local equilibrium state and monitor the time evolution of 
a* = cl/vq and P* y = P xy /p- Here, u = p/r] oc T 1 / 2 is an effective collision frequency, i] 
is the elastic shear viscosity. The simulations show that, after a transient regime, the ratio 
—P* y /a* reaches a constant value (independent of the shear rate), which can be identified 
as the shear viscosity in the linear hydrodynamic regime. Details of the simulation for dense 
fluids will be published elsewhere |25|] , and here only compare the Monte Carlo simulations 
for the shear viscosity with the Sonine approximations given by Eqs. (|28| ) and (|46|). 

In our simulations we have typically taken 10 5 particles and have averaged over 5 replicas. 
Since the thermal velocity Vq is not constant in the transient regime, we have taken a time- 
dependent time step given by 0M£/v o (t), where t = (v^vma 2 )" 1 is the mean free path. 
The simulation results are shown in Fig. |I[ In general, we observe a very good agreement 
between the predictions of the Chapman- Enskog theory in the first Sonine approximation 
and the simulation data. This agreement is similar to the one previously found in the free 



cooling case P0f . At a quantitative level, we see that the discrepancies between theory and 
simulation tend to increase as the dissipation increases, although these differences are quite 
small (less than 3%). As a final conclusion, it is important to remark that the agreement 
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extends over a wide range of values of the restitution coefficient (a > 0.6), indicating the 
reliability of the Sonine approximation for describing granular flows beyond the quasielastic 
domain. 



V. SUMMARY AND DISCUSSION 

In this paper we have addressed the derivation of the hydrodynamic equations of a 
granular gas from the Boltzmann kinetic theory. The system is heated by the action of 
"thermostatting" external forces which exactly compensate for cooling effects associated 
with the inelasticity of collisions. Two different types of thermostats have been considered: 

(a) an "anti-drag" force proportional to the particle peculiar velocity (Gaussian force), and 

(b) a stochastic force, which gives frequent kicks to each particle between collisions. The in- 
troduction of these thermostats has the advantage of avoiding the intrinsic time dependence 
of the homogeneous cooling state (unforced gas), but at the price of introducing unknown 
new effects induced by the external forcing. While in the homogeneous problem the results 
obtained with and without a Gaussian thermostat are completely equivalent (when one 
scales the particle velocity with the thermal velocity), this equivalence does not hold in the 
stochastic case and the non-Maxwellian properties of the distribution function (cumulants) 
are different from those obtained in the unforced case. Here, our goal has been to assess the 
influence of thermostats on the transport properties of the gas. 

The transport processes considered are those for a fluid with small spatial gradients of the 
hydrodynamic fields. For this reason, the Boltzmann equation has been solved perturbatively 
using an adaptation of the Chapman- Enskog method recently proposed for inelastic collisions 
PfJ]. By using similar procedures as those made in the free cooling case ||, we calculate the 
distribution function to first order in the gradients. Its use in the functionals for the pressure 
tensor and heat flux provides a representation of these as linear combinations of the gradients. 
The corresponding coefficients in these expressions are the shear viscosity r\ (defined in Eq. 
(p6|)), the thermal conductivity k and the new coefficient /i (both coefficients defined in 
Eq. (|27D). These transport coefficients are in general functions of the restitution coefficient 
a. Their expressions are given by Eqs. (|28|)-(|30"|) in the case of the Gaussian thermostat 
and Eqs. (f46"l) -(f48D in the case of the stochastic thermostat. A practical evaluation of these 
coefficients is possible by means of a Sonine polynomial approximation and the derivation 
and approximate results are not limited to weak inelasticity. 

The dependence of rj, k, and p on a has been illustrated in the case of hard spheres 
(d = 3). As Figs. 1]-^] show, the thermostats affect the transport properties since the 
discrepancies between the driven and free cooling results are quite significant. Although 
not widely recognized, the above conclusion illustrates the fact that generally the inclusion 
of an external force depending on the state of the system changes the apparent transport 
coefficients. This has been demonstrated as well for uniform shear flow with elastic collisions 
where thermostats are used to produce a steady state. In that case, the Navier-Stokes shear 
viscosity is unchanged, but nonlinear rheological properties are affected by the thermostat 
26fl . In the context of granular fluids, the effects already occur at the Navier-Stokes order. 



Notice that the above conclusion only affects to this type of external forcing mechanisms 
(termostats), since driving the system by shaking, vibration, and even the action of a weak 
external field (such as the gravity field) does not modify the transport coefficients of the 
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gas. Concerning the influence of dissipation on transport, we observe that in general the 
deviation from the functional form for elastic collisions is more significant in the unforced 
gas than in the driven cases. In particular, the coefficient fi (which is zero in the elastic 
limit) is clearly different from zero in the unforced gas (for instance, //* ~ 0.27 at a = 0.8) 
while is negligible in the Gaussian and stochastic cases (for instance, at a = 0.5, fi* ~ 0.064 
for the Gaussian force while /i* ~ 0.012 for the stochastic force). 

To check the accuracy of the Sonine approximation, we have numerically solved the 



Boltzmann equation by means of the Direct Simulation Monte Carlo (DSMC) method (19 
for a granular gas under uniform shear flow. To control inelastic cooling, a thermostat 
force is introduced in the system. In these conditions, the (apparent) Navier-Stokes shear 
viscosity can be measured directly in the long time limit just as for the case of elastic 
collisions. This simulation method has been recently proved to be an efficient way of 
measuring the shear viscosity of a moderately dense granular gas. The comparison carried 
out here in the low-density regime shows that the Chapman-Enskog results in the first Sonine 
approximation exhibit an excellent agreement with the simulation data, even for moderate 
dissipation (say, for instance a = 0.7). These results indicate clearly the reliability of the 
quantitative predictions for transport coefficients from the Chapman-Enskog method with 
small spatial gradients but including strong dissipation. 

As said before, a study of long-range correlations in a granular system fluidized by the 



random stochastic force (|13[ ) has been recently made ||10|| . In order to analyze the decay of 
fluctuations, it was assumed that the transport coefficients are equal to the corresponding 
quantities given by the elastic Enskog theory [|TJ. According to the results obtained in this 
paper for a heated granular gas in the low-density regime, although this assumption can 
be considered as a quite good approximation for the shear viscosity rj and the coefficient 
/i, this is not true for the thermal conductivity since k* clearly differs from 1, even for 
small inelasticity. In this context, it would be interesting to reexamine the conclusions 
obtained in Ref. JTUJ when the true Enskog transport coefficients are considered. As said 
in the Introduction, the Enskog equation for a gas of inelastic hard spheres (d = 3) in the 
absence of thermostats has been recently solved up to Navier-Stokes order. Taking into 
account these results Q, the extension of the calculations carried out in this paper to higher 
densities can be easily made. The final expressions of the corresponding Enskog transport 
coefficients in the three dimensional case are displayed in Appendix [B|. We plan to extend 
these expressions for a d dimensional granular fluid in the near future. 
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APPENDIX A: CHAPMAN-ENSKOG EXPANSION 

In the case of the Gaussian thermostat, the velocity distribution /W obeys the kinetic 
equation 
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{ c + l^W • v ) + l^W ' ( v/(0) ) = - + v • v ) /(0) ' (A1) 

with Df = df + u ■ V, £ is the linear operator given by 

£/(!) = _ + J[/W, /W]), (A2) 

and £W is a linear functional of J 1 - 1 ) defined as 

C (1) = --^ / dvmV* (J[/W, /«] + J[/(D, /(o)]) . (A3) 

The macroscopic balance equations to first order in the gradients are 

IT 

Dl 1] n = -nV ■ u, Df'ui = -(mn) _1 V^, D t (1) T = — — V ■ u. (A4) 



Use of these in Eq. (|A1|) yields 

[ C + ' V ) /W + l^Wv ' ( V/(0) ) = A • VlnT + B ■ Vlnn + QjV^j, (A5) 

where the expressions of A, B, and are the same as those obtained in the Appendix A 
of Ref. |H, i.e., 

A(V) - iv A . ( V /<«») - 1 A/<°>, (A6) 
B(V) = -V/<») - A/'", (A7) 

^CV) = A( Vj/m )_^.(v/m). (A8) 

The fact that Cy is traceless implies that the scalar £W = by symmetry. This is special to 
the low density Boltzmann equation, since at higher densities |3J there is a contribution to 
fi 1 ) proportional to V • u leading to a nonzero value of C^. This can be seen in Appendix 
[B|. Comparison with the kinetic equation obeying in the unforced case (Eq. (A5) of Ref. 
0) shows that both kinetic equations only differ in the terms d\ ; /i (unforced description) 
and (C (0) /2)(<9/<9V) ■ V/i (1) (driven case). These terms give different contributions to the 
transport coefficients. The solution to Eq. (|A5|) is of the form 



/(!) = A-V\nT + B-V\nn + C ij V i u j . (A9) 

Substituting Eq. (|A9|) into Eq. (|A1| ) and identifying coefficients of independent gradients, 
one gets the following three integral equations determining the unknowns A(V), B(V), and 
Cy(V): 
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( A 




B 


Cij J 





(A10) 



The expressions (p8|), (p9|) , and (^) for the transport coefficients follow directly from 
the integral equations ( |A1(J| ). The transport coefficients are defined as [[3] 



J dVZ^(V)C y (V), 



(d- l)(d + 2) 
/*=-^/^CV).B(V). 



(All) 



(A12) 



(A13) 



The thermal conductivity k and the coefficient /i can be easily obtained by multiplying the 
two first equations of ( |A10|) by S(V) and integrating over the velocity. The result is 



(*« - |c (0) ) « 



d + 2nT 



l + c) 



(A14) 



- |C {0) ) /i 



d + 2T 2 



-c. 



4 m 

Upon deriving these expressions, use has been made of the relation 



2 9V V J 



-Td T f {0) 



(A15) 



(A16) 



which follows from the temperature dependence of f(°\ The shear viscosity 77 can be obtained 
in a similar way by multiplying the third equation of ( |A10| ) by Dy-(V): 



(v v - C (0) ) 77 = nT. 



(A17) 



Equations (|A14|) , (|A15|) , and (|A17|) lead directly to expressions (|28| ) - (|30|) appearing in the 
text. 

The analysis in the case of the stochastic thermostat (|i~4|) is similar to that made above 
for the Gaussian one. To first order in gradients, one has 



2 m [dV 



/W = A- VlnT + B- Vlnn + C^V^. 



(A18) 



Here, the velocity dependence on the right side of Eq. ( |A18| ) is given by Eqs. ([A6|)-([A8|) and 
we have taken into account that £W = 0. The solution to Eq. (|A18|) is of the form ( |A9|) , 
where now the corresponding integral equations are 
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2 m [dV 



2 






( B 






{ di ) 













(A19) 



Now, we multiply Eq. (|A19[) by S(V) and Dij(V) and integrate over the velocity to get the 
expressions of k, /i, and 77: 



d + 2nT. 

= — (1 + c), 

2 m 



d + 2T 2 

c, 

4 m 

u v r] = nT. 



Here, we have made use of the results 

J dw S(V) • 



'jr 



A = 0, 



(A20) 

(A21) 
(A22) 

(A23) 



9V 



B = 0, 



(A24) 



/ dvA 3 (V) 

which follow directly from the solubility conditions (|I~9"D and 



(A25) 



APPENDIX B: RESULTS FOR A DENSE FLUID OF HARD SPHERES 

In this Appendix, we display the results derived for inelastic hard spheres (d = 3) in the 
framework of the Enskog equation when the gas is heated by the Gaussian and the stochastic 
thermostats. To first order in the gradients, the momentum and heat fluxes are 



iP = -V [ V<«i + VjUi - -6ijV • u ) - j5ij V ■ u, 



q 



(i) 



-kVT — /iVn, 



(Bl) 
(B2) 



where 7 is the bulk viscosity coefficient which vanishes in the low density limit. In a compact 
form, the transport coefficients can be written as 



7] — Tj 



1 + 



2t™*x(1 + a) 
15 



3 



(B3) 
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7 = ^VirmTn*nx(T{l + a) (l - ^ 



(B4) 



K = K 



1 + 



7m*x(l + a) 



+ 3V^ Xn ° l 1 + ^J 



(B5) 



/X = fi 



1 + 



7rn*x(l + a) 



(B6) 



where x is the pair correlation function at contact, n* = na 3 is a reduced density, and the 
superscript k denotes the contributions from the kinetic parts of the fluxes [|J] . These kinetic 
parts are given by 



V k = nT (^-6,C (0) ) 1 



1 — — (1 + a)(l — 3a)nn*Y 
15 



(B7) 



2 m v y 



1 + c ( H — 7m x 



+ — 7m*x(l + a) 2 



2a - 1+ 2( 1+ «)-^y 



(B8) 



i f5nT 



1 + a 



1 + n<9 n In ( 1 H i — ™*X 



nT / 1 



-Trn x - „ ... 
«(1 — a) + - fa(l — a) + - ) ' 



(B9) 



In these equations, = —1, and b K = 6 M = — | in the case of the Gaussian thermostat while 



b v — b K — b/j, — for the stochastic thermostat. 



Up to first order in gradients, the cooling rate ( is given by 



--7m* v + — z/ f 1 + — ) y/3 
3 A 32 ° V 64 M 



l-a 2 )V-u, (BIO) 



where the expression of the quantity (3 depends on the thermostat used. In the Gaussian 
case, one has 



(2/45) A - (c/6)(l + a) (5 - 3a) 
I / 7 -2C(°)-(5c/64)(1-q2) (l + |)x 



(Bll) 



while in the stochastic case one gets 



(3 = 7in*x 



(2/45)A-(c/3)(l + a) 



/a 



(B12) 



Here, we have introduced the quantities 
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A = -(l- a 2 ) 5a 2 + 4a-l + 
8 V ; 1 12 



c 159a + 3a 2 -19a- 15a 3 ' 



1 — a 



(B13) 



1 + a 
^8" 



128 - 96a + 15a 2 - 15a 3 + — fl5a 3 - 15a 2 + 498a - 434) 

64 V > 



(B14) 



To get the explicit dependence of the transport coefficients and the cooling rate on the 
reduced density n*, one can take for instance the Carnahan-Starling approximation for \ 
given by 



X 



12 



7rrT 



1 — ^itn* 



(B15) 



It is easy to check that all results of Appendix reduce to those presented in the text for 
the low-density limit (n* = 0). 
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FIGURES 



FIG. 1. Reduced shear viscosity 77* = 77/770 as a function of the restitution coefficient a for 
a three-dimensional system. The solid lines refer to the theoretical expressions derived in the 
unforced case (a), in the Gaussian case (b), and in the stochastic case (c). The symbols are the 
results obtained from the Direct Simulation Monte Carlo method in the Gaussian (circles) and 
stochastic (triangles) cases. 

FIG. 2. Reduced thermal conductivity k* = k/ko as a function of the restitution coefficient 
a for a three-dimensional system in the unforced case (a), in the Gaussian case (b), and in the 
stochastic case (c). 

FIG. 3. Reduced coefficient fi* = nfi/Tno as a function of the restitution coefficient a for a 
three-dimensional system in the unforced case (a), in the Gaussian case (b), and in the stochastic 
case (c). 
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